# Code to overlay the 3 graphs
# @author Sahir Shahryar

n50 <- read.csv(file = "POLSData/bayes.50n.csv")
n150 <- read.csv(file = "POLSData/bayes.150n.csv")
n250 <- read.csv(file = "POLSData/bayes.250n.csv")

iterations <- c(1:20, 30, 40, 50)

cols <- c("Beta 0", "Beta 1", "Beta 2", "Beta 3", "Beta 4", "Beta 5", "Partial Sill", "Range", "Nugget/Partial Sill")
names <- c("Beta0", "Beta1", "Beta2", "Beta3", "Beta4", "Beta5", "PartialSill", "Range", "RelativeNugget")

for (i in 1:length(cols)) {
	col <- cols[i]
	max <- max(n50[, i + 1], n150[, i + 1], n250[, i + 1])

	jpeg(paste0("bayes", names[i], "combined.jpg"))
	par(las = 1)
	plot(x = iterations, y = n50[, i + 1], col = "#FF0000", ylim = range(c(0, max * 1.3)), xlim = range(c(1, 50)), xlab = "Iterations", ylab = "Percent Mean Absolute Difference", main = col)
	points(x = iterations, y = n50[, i + 1], col = "#FF0000", pch = 15)
	points(x = iterations, y = n150[, i + 1], col = "#00FF00", pch = 17)
	points(x = iterations, y = n250[, i + 1], col = "#0000FF", pch = 20)
	lines(x = iterations, y = n50[, i + 1], col = "#FF0000", lty = 6)
	lines(x = iterations, y = n150[, i + 1], col = "#00FF00", lty = 5)
	lines(x = iterations, y = n250[, i + 1], col = "#0000FF", lty = 1)

	legend(1, max * 1.3, legend = c("n = 50", "n = 150", "n = 250"), col = c("#FF0000", "#00FF00", "#0000FF"), lty = c(6, 5, 1), pch = c(15, 17, 20))

	dev.off()
}
